#ifndef _CGCM_DATATYPES_
#define _CGCM_DATATYPES_

#include <vector>
#include <string>
#include <stack>
#include <assert.h>
#include <list>
#include "precision.h"
#include "qgcmmath.h"
#include <complex>

using namespace std;

/*****************************************************************************************************
  General description:

  Declares structs and classes used by model's compiler, p-system (ODE solver), 
  and s-system (PDE solver)
 ******************************************************************************************************/

#ifdef DEBUG // for use with gdb (output data to the remote output terminal: rout)
FILE * getsetrout(const char *);
FILE * getrout();
void rclear();
void rprint(const char *);
#endif


// ***  --> typedefinitions
typedef char INDEX;
typedef PRECISIONTYPE REAL;
typedef complex<REAL> qgcomplex; // short hand notation


//***************************************************************************************************/
// rank class: indexing object needed for multi-dimensional arrays
//***************************************************************************************************/
template <class T, int s> 
struct rank {
	T list[s];
	bool operator==(const rank&) const;
	T & operator[](size_t);
	const T & operator[](size_t) const;
	void offset0();
	void operator=(const rank&);
	void operator--();
	// void operator=(const double &);
	bool operator<(const rank&) const;
	bool operator>(const rank&) const;
	rank();
	rank(const T , const T );
	rank(const T , const T , const T );
	rank(const T , const T , const T , const T );
	rank(const T , const T , const T , const T, const T );
	rank(const T , const T , const T , const T , const T , const T);
	rank(const T , const T , const T , const T , const T, const T, const T, const T, const T, const T, const T );
	rank(const T , const T , const T , const T , const T, const rank<T, 3> &, const rank<T, 3> &);
	rank(const rank &);
#ifdef DEBUG
	string print() const;
#endif
};

// force instantiation of the following templates
//-------------------------------------------------------------------------------------------
typedef rank<INDEX,2> rank2;
typedef rank<INDEX,3> rank3;
typedef rank<INDEX,4> rank4;
typedef rank<INDEX,5> rank5;
typedef rank<INDEX,6> rank6;
typedef rank<INDEX,11> rank11;


template <class T, int s>
class rankIter {
	private:
		rank<T,s> bounds;
		rank<T,s> value;
	public:
		bool operator==(const rankIter&) const;
		bool operator!=(const rankIter&) const;
		rankIter & operator ++();
		rank<T,s> & operator*();
		bool operator<(const rank<T,s> &) const;
		bool operator>(const rank<T,s> &) const;
		rankIter(const rankIter &);
		rankIter(const rank<T,s> & bounds);
		rankIter(const rank<T,s> & bounds, const rank<T,s> & initialValue);
};
// force instantiation of the following templates
//-------------------------------------------------------------------------------------------
typedef rankIter<INDEX,6> rank6iter;
typedef rankIter<INDEX,4> rank4iter;

//***************************************************************************************************/
// Multi-dimensional array class: Basically a map between a 1d array (which is how the data is stored)
//   and a nth dimensional array object, which is how the physics is often represented.
//***************************************************************************************************/
template <class T, int s>
class multidimarray { // number of dimensions of array is fixed by second parameter
	T * data;
	rank<INDEX, s> edges;
	public:
	multidimarray();
	multidimarray(const multidimarray &);
	multidimarray & operator=(const multidimarray &);
	multidimarray(const T &, const rank<INDEX, s> &  );
	~multidimarray();
	void resize(const rank<INDEX, s> &, const T &);
	T & operator[](size_t);
	const T & operator[](size_t) const;
	T & operator[](const rank<INDEX, s> &); 
	const T & operator[](const rank<INDEX, s> &) const;
	const rank<INDEX,s> & getedges() const;
	size_t getlinearsize() const;
	inline const T * getptr() const;
	T * getdata(); 
};


//***************************************************************************************************/
// regset: data blocks used to store model data
//         There are two: 1) primary storage, 2) intermediate storage used for temporarily storing
//                           provisional values during forward time integration for corrector
//			     schemes.
//***************************************************************************************************/
typedef multidimarray< qgcomplex ,4> array4d;
struct regset {
	mutable qgcomplex * data; // data always points to either the primary or intermediate block
	// its value can be changed without changing the values of the data itself
	// so the element is mutable
	qgcomplex * primary_data;
	qgcomplex * intermediate_data;
	rank4 edges;
	int lsize;

	regset();
	regset & operator=(const regset  &);
	regset(const regset &);
	bool allocate(const rank4 & edges);
	void deallocate();
	void setintermediate() const;
	void setprimary() const;
	void zeromemory();
	qgcomplex & operator[](unsigned int index); // index data as 1d array
	const qgcomplex & operator[](unsigned int index) const; // index data as 1d array
	const qgcomplex & operator[](const rank4 & index) const; // index by rank object
	qgcomplex & operator[](const rank4 & index); // index by rank object
	qgcomplex * get_data();
	const qgcomplex * get_data() const;
	qgcomplex * get_primary_data();
	qgcomplex * get_intermediate_data();
	int size();
	const int size() const;
#ifdef DEBUG
	void dump() const;
	void dumpf(const char * fname) const;
	void dumph(FILE *) const; 
#endif
	~regset();
};
//****************************************************************************
class twotimevar {
	const rank4 edges; // edges on all the regsets
	vector< regset > elements; // forming the triangle
	int maxt; 
	private: // prohibit compiler from supplying an implicit contructor:
	twotimevar(); // no empty constructor allowed

	public:
	twotimevar(regset & initdata);
	void addtime();
	void addtimes(const vector<int> & list);
	void addeqtime();

	bool does_element_exist( unsigned int, unsigned int ) const;
	regset & getvar(unsigned int, unsigned int);
	const regset & getvar(unsigned int, unsigned int) const;
	const regset & getvar_withcgcmexception(unsigned int, unsigned int) const;
	void delvar(unsigned int, unsigned int);
	const int getmaxtime() const; 
	// allow the user to eliminate integration history if they desire it
	// this could be dangerous in fully nonlinear runs
};

//***************************************************************************************************/
// parameter: stores  model parameters
//***************************************************************************************************/
class parameter {
	public:
		enum type { Real, String, Integer };
	private:
		union {
			vector<string> * svalues;
			vector<int> *    ivalues;
			vector<REAL> *   fvalues;
		} values;
		string name;
		type ptype;
		bool constant;
	public:
		const vector<string> & getsData() const;
		const vector<int> & getiData() const;
		const vector<REAL> & getData() const;
		const type getType() const;
		void setValue(REAL);
		void setValue(REAL, int);
		void setValue(int);
		void setValue(int, int);
		void setValue(const char *);
		void setValue(const char *, int);
		const REAL & getfValue() const;
		const int  & getiValue() const;
		const char * getsValue() const;
		const char * GetName() const;
		const size_t GetSize() const;
		void assignData( const vector<REAL> & );
		void assignData( const vector<int> & );
		void assignData( const vector<string> & );
		void MakeConstant();
		bool IsConstant() const;

		void resize( int );

		parameter(); // default constructor ....
		parameter( const string &, const int & );
		parameter( const string &, const REAL & );
		parameter( const string &, const string & );
		parameter( const string &, const vector<int> & );
		parameter( const string &, const vector<REAL> & );
		~parameter();
};
//***************************************************************************************************/
// local exception class, might require further complexity in the future?
class cgcmexception : public exception {
	public:

		cgcmexception() throw();
		cgcmexception(const char *, int) throw ();
		cgcmexception(const string &, int) throw();
		~cgcmexception() throw();

		const char * what() const throw();
		int id() const throw();
	private:
		const int rvalue;
		const string reason;
};

#endif
